%----------------------------
% clear plots
 close all;
%---------------------------

// Christoph Gortz and John Tsoukalas
// News and Financial Intermediation in Aggregate Fluctuations, Review of Economics and Statistics

// Two sector model with price and wage rigidities with permanent and transitory 
// technology shocks and separate specification of Multi Factor Productivity in I anc C sectors
// Capital is immobile across sectors and there are two investment and capital utilization decisions
// intratemporal investment adjustment cocts are incorporated a-la Huffman and Whynne
// This version does not include an Investment-specific shock v_l and a TFP transitory shock to cons sector z_l!
// This model includes a banking sector a la Gertler and Karadi (2010) 

var   labi labc lab kc ki kpi kpc mcc mci zcapi zcapc rk pi lam phifi phifc c invei invec inve y pinfc pinfi w r 
z v b g qs ms sw ewma z_l v_l dy dc dinve    dpi   dw labobs robs pinfobsc pinfobsi spreadobsc spreadobsi rrate //spread dcredit
epinfmac epinfmai spinfc spinfi spi spkc spki snc sni
 rbc rbi  qc qi epkc1aux epkc2aux epkc3aux epkc4aux epkc8aux epki1aux epki2aux epki3aux epki4aux epki8aux
ez1aux ez4aux ez8aux ev1aux ev4aux ev8aux tfp etfp4aux etfp8aux dtfp ezl4aux ezl8aux evl4aux evl8aux me_tfpc 
cthetab 
;


varexo ez ev eb em epinfc epinfi ew enc eni eg epkc epki eqs evl ezl espi epkc1 epkc2 epkc3 epkc4 epkc8 epki1 epki2 epki3 epki4 epki8 
ez1 ez4 ez8 ev1 ev4 ev8 etfp etfp4 etfp8 ezl4 ezl8 evl4 evl8 emetfpc emetfpi ethetab;
 
//esp1, esp2 are measurement errors

parameters constepinfc constepinfi constebeta cmaw cmapi cmapc  alfac alfai
czcapi czcapc czcap csadjcost cdeltai cdeltac csigma chabb  
cindw cprobw cindpc cprobpc cindpi cprobpi csigl cfc
crpi crdy cry crr 
crhoz crhov crhob crhog crhoqs crhoms crhow crhozl crhovl crhospkc crhospki
cminrho
cg cga 
cgv clandapaux clandawaux
Lss nvalueic 
crhopinfc crhopinfi crhospi
cthetabss crhothetab clamb crhosnc crhosni ctau cvarthetac cvarthetai psicss psiiss lridata lrcdata //equityss//cvaromega
grconstepinfc crhotfp crhotfp_mec invshare

// Steady state parameters

X1aux crk cpi cw KcLc KiLi Lc Li css iss Kc Ki K Kibar Kcbar yss iiss icss
// Parameters needed to calculate ss parameters
cbeta ga 
gv clandap clandaw
// Steady state parameters resulting from financial intermediaries
z1css z1iss z2css z2iss rbcss rbiss rss lrcss lriss nucss nuiss qcss qiss sncss sniss scss siss;


// fixed parameters

cdeltai =.025;   // depreciation rate I sector capital
cdeltac =.025;   // depreciation rate C sector capital
alfac=.3;
alfai=.3;

// SS restriction is cbeta = pinfobs/R, for US data (1990Q2:2011Q1) this gives cbeta = 0.9974
// from this you can get implied value for constebeta = (cbeta^-1 -1)*100
constebeta = 0.2607;   // equivalent to cbeta = 0.9974, constistent with 1990Q2-2011Q1 US data
constepinfc = 0.6722;  // constistent with 1990Q2-2011Q1 US data
grconstepinfc = constepinfc/100 + 1;
constepinfi = 0.0245;  // constistent with 1990Q2-2011Q1 US data
cg = 0.177;            // ss nominal government spending to nominal output ratio, constistent with 1982Q4-2011Q1 US data
nvalueic = 0.399;      // average value of nominal investment to consumption ratio, constistent with 1982Q4-2011Q1 US data
Lss = 1;               // ss of total hours (by choosing appropriate value of the weight of leisure in utility)
cmapc =    0;          // MA price mark-up
cmapi =    0;          // MA price mark-up
cmaw  =   0;           // MA wage mark-up
cminrho = 0.0000001;//0.3578;      // intrtemporal inv. adjustment cost, zero equals crho = -1
csigma=1.0;
cfc=1.5;




// estimated parameters initialisation

cga = 0.141; 
cgv = 0.434;
csadjcostc= 4;	// investment adj. cost parameter
csadjcosti= 4;	// investment adj. cost parameter

chabb=    0.5;      // consumption habit 
cprobw=   0.01;     // wage calvo contract
cindw=    0.01;      // wage indexation
csigl=    2;      // inverse labour supply ela
cprobpc=   0.01;    // consumption sector price calvo probability 
cindpc=    0.01;     // consumption sector price indexation
cprobpi=   0.01;    // investment sector price calvo probability 
cindpi=    0.01;     // investmnt sector price indexation
czcapi=    5;    // variable I-capital utilization 
czcapc=    5;    // variable C-capital utilization 
czcap =    5.75;
crpi=     1.7;      // taylor rule inflation 
crr=      0.6;        // taylor rule inertia  
cry=      0.0;        // taylor rule output gap 
crdy=     0.25;        // taylor rule output gap growth
clandawaux = 1.01-1; // ss wage mark up
clandapaux = 1.01-1; // ss price mark up

crhoz=    0.4;      // consumption sector productivity AR(1) --- growth shock
crhov=    0.4;        // investment sector productivity AR(1) --- growth shock
crhozl=    0.0;     // consumption sector productivity AR(1) --- stationary shock
crhovl=    0.0;     // investment sector productivity AR(1) --- stationary shock
crhob=    0.6;      // preference shock AR(1)
crhog=    0.6;      // gov spending shock AR(1)
crhoqs=   0.0;      // installation shock AR(1)
crhoms=   0.0;      // monetary policy shock AR(1)
crhopinfi= 0.6;     // price mark-up shock AR(1) in I sector
crhopinfc= 0.6;     // price mark-up shock AR(1) in C sector
crhow=    0.6;      // wage mark-up shock AR(1)
crhospkc = 0.6;     // quality of physical capital shock in C sector
crhospki = 0.6;     // quality of physical capital shock in I sector
crhospi = 0.0;      // persistence of inflation objective shock
crhosnc = 0.0;      // Persistence of shock to intermediarys equity capital in C-sector
crhosni = 0.0;      // Persistence of shock to intermediarys equity capital in I-sector
crhothetab = 0;     // Persistence of equity shock
crhotfp = 0.0;      // Persistence common TFP
crhotfp_mec = 0;    // Persistence measurement error TFP

// Fixed parameters specific to financial intermediaries
cthetabss = 0.96;      // Fraction of bankers that survive at following period (in SS)
lridata = 5.47;         // Leverage ratio implied by the data for I sector
lrcdata = lridata;      // Leverage ratio implied by the data for C sector
cvarthetac = 0;         // Response of gov's credit to external spread in C sector
cvarthetai = cvarthetac;// Response of gov's credit to external spread in I sector
psicss = 0;             // Steady state fraction assets intermed by CB in C sector
psiiss = psicss;        // Steady state fraction assets intermed by CB in I sector
ctau = 0.0;             // Government's cost of providing profit
//equityss = 0.5;         // C-sector equity over total equity in steady state
// Central bank as an intermediary is shut down for psicss = psiss = ctau = cvarthetac = cvarthetai = 0


// Steady state parameters
// Initial values according to steady state with rho = -1
// parameters needed for this section:
cbeta = 1/(1+constebeta/100);
ga = cga/100 ;
gv = cgv/100 ;
clandap = 1 + clandapaux;
clandaw = 1 + clandawaux;
//constepinfi  = constepinfc + (cga+ (alfac-1)/(1-alfai)*cgv);

invshare = 0.24;



model (linear); 

// Auxiliary parameter
//# crho = -(cminrho + 1);  // rho = -rho
# crho = 1/(cminrho - 1);   // rho = -rho
//# cga = constepinfi - constepinfc -  (alfac-1)/(1-alfai)*cgv;
//# ga = cga/100;




       

// sticky price - wage economy



          // consumption and investment sector output and cost minimization
          
         inve = clandap*(alfai*ki + (1-alfai)*labi + v_l + tfp);
          inve = (iiss^-crho + icss^-crho)^-1*iiss^-crho*invei + (iiss^-crho + icss^-crho)^-1*icss^-crho*invec;
          
          //inve = (iiss/iss)*invei + (icss/iss)*invec;
          c  + crk*Kibar/css*exp(-(1/(1-alfai))*gv)*zcapi + crk*Kcbar/css*exp(-(1/(1-alfai))*gv)*zcapc = clandap*(alfac*kc + (1-alfac)*labc + z_l + tfp);

	      mcc =  alfac*rk+(1-alfac)*w - z_l - tfp;
          mci =  alfai*rk+(1-alfai)*w - pi - v_l - tfp;
          rk - w = labc - kc;
          rk - w = labi - ki;

         // consumption and investment sector optimal pricing --- Philips curves
          
           pinfc -spi =  (cbeta/(1+cindpc*cbeta))*(pinfc(+1)-spi) + (cindpc/(1+cindpc*cbeta))*(pinfc(-1)-spi) 
               +( (1-cprobpc)*(1-cbeta*cprobpc)/((1+cindpc*cbeta)*cprobpc) )* mcc  +  spinfc ;
               
           pinfi -spi =  (cbeta/(1+cindpi*cbeta))*(pinfi(+1)-spi) + (cindpi/(1+cindpi*cbeta))*(pinfi(-1)-spi) 
               +( (1-cprobpi)*(1-cbeta*cprobpi)/((1+cindpi*cbeta)*cprobpi) )* mci  + spinfi ; 

          
          // normalization of shocks:
          // spinf* = ( (1-cprobpi)*(1-cbeta*cprobpi)/((1+cindpi*cbeta)*cprobpi) )* spinf ; 
            
        // Households
        
        lam = (chabb*cbeta*exp(ga+(alfac/(1-alfai))*gv)/((exp(ga+(alfac/(1-alfai))*gv)-chabb*cbeta)*(exp(ga+(alfac/(1-alfai))*gv)-chabb)))*c(+1)
        - ((exp(2*(ga+(alfac/(1-alfai))*gv)) + chabb^2*cbeta)/((exp(ga+(alfac/(1-alfai))*gv)-chabb*cbeta)*(exp(ga+(alfac/(1-alfai))*gv)-chabb)))*c
        + (chabb*exp(ga+(alfac/(1-alfai))*gv)/((exp(ga+(alfac/(1-alfai))*gv)-chabb*cbeta)*(exp(ga+(alfac/(1-alfai))*gv)-chabb)))*c(-1)
        + ( exp(ga+(alfac/(1-alfai))*gv)*chabb*cbeta*crhoz - chabb*exp(ga+(alfac/(1-alfai))*gv) )/((exp(ga+(alfac/(1-alfai))*gv)-chabb*cbeta)*(exp(ga+(alfac/(1-alfai))*gv)-chabb))*z
	    + b
        +  (alfac/(1-alfai))*( exp(ga+(alfac/(1-alfai))*gv)*chabb*cbeta*crhov - chabb*exp(ga+(alfac/(1-alfai))*gv) )/((exp(ga+(alfac/(1-alfai))*gv)-chabb*cbeta)*(exp(ga+(alfac/(1-alfai))*gv)-chabb))*v;
        
	    //normalization of shocks:
        //b* =  ( exp(ga+(alfac/(1-alfai))*gv) - chabb*cbeta*crhob )/(exp(ga+(alfac/(1-alfai))*gv)-chabb*cbeta)*b

        lam = r + lam(+1) - z(+1) - v(+1)*alfac/(1-alfai) - pinfc(+1);
        
        phifi = (1-cdeltai)*cbeta*exp(-(1/(1-alfai))*gv)*( phifi(+1) - 1/(1-alfai)*v(+1) + spki(+1) )
        + (1 - (1-cdeltai)*cbeta*exp(-(1/(1-alfai))*gv))*( lam(+1) - 1/(1-alfai)*v(+1) + rk(+1) + zcapi(+1) + spki(+1) );

        phifc = (1-cdeltac)*cbeta*exp(-(1/(1-alfai))*gv)*( phifc(+1) - 1/(1-alfai)*v(+1) + spkc(+1))
        + (1 - (1-cdeltac)*cbeta*exp(-(1/(1-alfai))*gv))*( lam(+1) - 1/(1-alfai)*v(+1) + rk(+1) + zcapc(+1) + spkc(+1));

        
        rk = czcapi*zcapi;
        rk = czcapc*zcapc;

        //zcapi =  (1/(czcap/(1-czcap)))* rki ;
        //zcapc =  (1/(czcap/(1-czcap)))* rkc ;
        
        lam = -pi + phifi + qs - exp(2*((1/(1-alfai))*gv))*csadjcost*(invei - invei(-1) + 1/(1-alfai)*v) 
        + cbeta*exp(2*((1/(1-alfai))*gv))*csadjcost*( invei(+1) - invei +  1/(1-alfai)*v(+1) ) -(1+crho)*((iiss^-crho + icss^-crho)^-1* (icss^-crho*invec+iiss^-crho*invei)-invei);

        lam = -pi + phifc + qs - exp(2*((1/(1-alfai))*gv))*csadjcost*(invec - invec(-1) + 1/(1-alfai)*v) 
        + cbeta*exp(2*((1/(1-alfai))*gv))*csadjcost*( invec(+1) - invec +  1/(1-alfai)*v(+1) ) -(1+crho)*((iiss^-crho + icss^-crho)^-1* (icss^-crho*invec+iiss^-crho*invei)-invec);
        

       
	    ki =  kpi(-1)+ zcapi + spki - 1/(1-alfai)*v;

        kc =  kpc(-1)+ zcapc + spkc - 1/(1-alfai)*v;
 
        //kpi =  (1-cdeltai)*exp(-(1/(1-alfai))*gv)*( kpi(-1) + spki - (1/(1-alfai))*v ) + ( 1 - (1-cdeltai)*exp(-(1/(1-alfai))*gv))*(qs + invei);

        //kpc =  (1-cdeltac)*exp(-(1/(1-alfai))*gv)*( kpc(-1) + spkc - (1/(1-alfai))*v ) + ( 1 - (1-cdeltac)*exp(-(1/(1-alfai))*gv))*(qs + invec);

        Kcbar*kpc + Kcbar*kpi = (1-cdeltac)*exp(-(1/(1-alfai))*gv) * (Kcbar*kpc(-1) + Kcbar*spkc - 1/(1-alfai)*Kcbar*v) + (1-cdeltai)*exp(-(1/(1-alfai))*gv) * (Kibar*kpi(-1) + Kibar*spki - 1/(1-alfai)*Kibar*v) + (Kcbar+Kibar - (1-cdeltac)*exp(-(1/(1-alfai))*gv)*Kcbar - (1-cdeltai)*exp(-(1/(1-alfai))*gv)*Kibar) * (qs + inve);

	      w =  (1/(1+cbeta))*w(-1)
               +(cbeta/(1+cbeta))*w(+1)
               - ( (1-cprobw*cbeta)*(1-cprobw)/(cprobw*(1+cbeta)*(1+csigl*(1+1/(clandawaux)))))* (w - csigl*lab - b + lam)
               + (cindw/(1+cbeta))*(pinfc(-1)-spi)
               - ((1+cbeta*cindw)/(1+cbeta))*(pinfc-spi)
               + (cbeta/(1+cbeta))*(pinfc(+1)-spi)
               + (cindw/(1+cbeta))*( z(-1) + alfac/(1-alfai)*v(-1) )
               - (1+cbeta*cindw-crhoz*cbeta)/(1+cbeta)*z
               - ((1+cbeta*cindw-crhov*cbeta)/(1+cbeta))*(alfac/(1-alfai))*v
               + sw ;
               
      // normalization of shocks:
      //       sw* = ( (1-cprobw*cbeta)*(1-cprobw)/(cprobw*(1+cbeta)*(1+csigl*(1+1/(clandaw-1)))))* sw 
               
      lab = Lc/Lss*labc + Li/Lss*labi;
               
      y = (css/yss)*c+(cpi*iss/yss)*(inve+ pi) + g ;
      //(1-cg)*y = (css/yss)*c+(cpi*iss/yss)*(inve+ pi) + (1-cg)*g + ctau/yss* ( psicss*qcss*scss*(psic+qc+sc) + psiiss*qiss*siss*(psii+qi+si) );         
        
      // Taylor rule         
	  //    r =  crpi*(1-crr)*pinfc
      //         +cry*(1-crr)*(y-yf)     
      //         +crdy*(y-yf-y(-1)+yf(-1))
      //         +crr*r(-1)
      //        +ms  ;

          r =  crpi*(1-crr)*pinfc
               +cry*(1-crr)*(pinfc - pinfc(-1))*0
               +crdy*(1-crr)*(y-y(-1))
               +crr*r(-1)
              +ms  ;

    //  Taylor rule similar to Smets Wouters 2005 to account for trend in interest rate observable
//        r =  spi + (1-crr)*( crpi*(pinfc(-1)-spi(-1)))
//               +(1-crr)*cry*(y-yf)     
//               +crdy*(y-yf-y(-1)+yf(-1))
//               +crr*(r(-1)-spi(-1))
//               +ms  ;
               
    pinfi - pinfc  = pi - pi(-1); 


   //dpi + (constepinfc - constepinfi) = pinfi  - pinfc;
// this implies the restriction
// constepinfi - constepinfc = (cga+ (alfac-1)/(1-alfai)*cgv)
// so initial parameter values need to respect this
// in order to solve the SS 

    
        // Stochastic return on assets for bank in C sector
    	rbc = 1/(crk+qcss*(1-cdeltac)) * ( crk*(rk+zcapc) + qcss*(1-cdeltac)*(qc) ) -qc(-1) +spkc + z - (1-alfac)/(1-alfai)*v;

        // Stochastic return on assets for bank in I sector
    	rbi = 1/(crk+qiss*(1-cdeltai)) * ( crk*(rk+zcapi) + qiss*(1-cdeltai)*(qi) ) -qi(-1) +spki + z - (1-alfac)/(1-alfai)*v;

        // Price of capital in C sector
		qc = exp(2*1/(1-alfai)*gv)*csadjcost * ( invec-invec(-1)+1/(1-alfai)*v) - cbeta*exp(2*1/(1-alfai)*gv)*csadjcost * (invec(+1)-invec+1/(1-alfai)*v(+1)) - qs + pi + (1+crho) * ( (iiss^(-crho)+icss^(-crho))^(-1) * (icss^(-crho)*invec + iiss^(-crho)*invei) - invec );

        // Price of capital in I sector
        qi = exp(2*1/(1-alfai)*gv)*csadjcost * ( invei-invei(-1)+1/(1-alfai)*v) - cbeta*exp(2*1/(1-alfai)*gv)*csadjcost * (invei(+1)-invei+1/(1-alfai)*v(+1)) - qs + pi + (1+crho) * ( (iiss^(-crho)+icss^(-crho))^(-1) * (icss^(-crho)*invec + iiss^(-crho)*invei) - invei );

       

        
 



    // Shock processes:
    // TFP growth shock to cons sector

        z = crhoz*z(-1)  + ez  + ez1aux(-1) + ez4aux(-4) + ez8aux(-8);
        ez1aux = ez1;
        ez4aux = ez4;
        ez8aux = ez8;
    
    // TFP stationary shock to cons sector
      
        z_l = crhozl*z_l(-1)  + ezl + ezl4aux(-4) + ezl8aux(-8);
        ezl4aux = ezl4;
        ezl8aux = ezl8;
         
   // Investment-specific sector growth shock

        v = crhov*v(-1)  + ev + ev1aux(-1) + ev4aux(-4) + ev8aux(-8);
        ev1aux = ev1;
        ev4aux = ev4;
        ev8aux = ev8;      

    // Investment-specific stationary shock 
      
         v_l = crhovl*v_l(-1)  + evl + evl4aux(-4) + evl8aux(-8);
        evl4aux = evl4;
        evl8aux = evl8;

   
   //  Preference, government, installation
   //  monetary policy, price markup, wage markup
	           
        b = crhob*b(-1) + eb;
	      
        g = crhog*g(-1) + eg;
        
        qs = crhoqs*qs(-1) + eqs;
	     
        ms = crhoms*ms(-1) + em;
	    
        spinfc = crhopinfc*spinfc(-1) + epinfmac - cmapc*epinfmac(-1);
	          epinfmac=epinfc;

        spinfi = crhopinfi*spinfi(-1) + epinfmai - cmapi*epinfmai(-1);
	          epinfmai=epinfi;
	    
        sw = crhow*sw(-1) + ewma - cmaw*ewma(-1) ;
	          ewma=ew; 
    // Shock to inflation objective
	      spi = crhospi*spi(-1) + espi;

    
    // Shock to bank's equity capital in C and I sector
       	snc = crhosnc*snc(-1) + enc;		
        sni = crhosni*sni(-1) + eni;
        //sn = crhosn*sn(-1) + en;	

    // Shock to the quality of physical capital in C and I sector
        spkc = crhospkc*spkc(-1) + epkc + epkc1aux(-1) + epkc2aux(-2) + epkc3aux(-3) + epkc4aux(-4) - epkc8aux(-8);
        epkc1aux = epkc1;
        epkc2aux = epkc2;
        epkc3aux = epkc3;
        epkc4aux = epkc4;
        epkc8aux = epkc8;
        spki = crhospki*spki(-1) + epki + epki1aux(-1) + epki2aux(-2) + epki3aux(-3) + epki4aux(-4) + epki8aux(-8);
        epki1aux = epki1;
        epki2aux = epki2;
        epki3aux = epki3;
        epki4aux = epki4;
        epki8aux = epki8;

// TFP shock common to both sectors
        tfp = crhotfp*tfp(-1) + etfp + etfp4aux(-4) + etfp8aux(-8);
        etfp4aux = etfp4;
        etfp8aux = etfp8;

// TFP measurement error shock  both sectors (with aggregate TFP as observable just use me_tfpc)
        me_tfpc = crhotfp_mec*me_tfpc(-1) + emetfpc ;

    // Shock to bank's equity
       cthetab = crhothetab*cthetab(-1) + ethetab;


// measurement equations
/*
dy=y-y(-1) + z + alfac/(1-alfai)*v + (cga + alfac/(1-alfai)*cgv) ;
dc=c-c(-1) + z + alfac/(1-alfai)*v + (cga + alfac/(1-alfai)*cgv)  ;
dinve=inve-inve(-1) + 1/(1-alfai)*v + ( 1/(1-alfai)*cgv) ;
dpi = pi - pi(-1) + z + (alfac-1)/(1-alfai)*v + (cga+ (alfac-1)/(1-alfai)*cgv) ;
dw=w-w(-1) + z + alfac/(1-alfai)*v + (cga + alfac/(1-alfai)*cgv)  ;
pinfobsc = 1*(pinfc) + constepinfc ;
pinfobsi = 1*(pinfi) + constepinfi ;
robs =    r + log(constepinfc/100 + 1) - log(cbeta) + (cga + alfac/(1-alfai)*cgv) ;
labobs = lab + log(Lss) ;
rrate = r - pinfc(+1); 
spreadobsi = rbi(+1)+log(rbiss)-robs;
spreadobsc = rbc(+1)+log(rbcss)-robs;
dequity = exp(ga + alfac/(1-alfai)*gv) * ( ( (qcss*scss/lrcss)/(qcss*scss/lrcss+qiss*siss/lriss) ) *(nc-nc(-1)) 
+ (1- ( (qcss*scss/lrcss)/(qcss*scss/lrcss+qiss*siss/lriss) ) )*(ni-ni(-1)) + z + alfac/(1-alfai)*v +cga + alfac/(1-alfai)*cgv);
*/
// NEW measurement equations (using demeaned data)

dy=y-y(-1) + z + alfac/(1-alfai)*v   ;
dc=c-c(-1) + z + alfac/(1-alfai)*v   ;
dinve=inve-inve(-1) + 1/(1-alfai)*v  ;
dpi = pi - pi(-1) + z + (alfac-1)/(1-alfai)*v  ;
dw=w-w(-1) + z + alfac/(1-alfai)*v   ;
pinfobsc = 1*(pinfc)  ;
pinfobsi = 1*(pinfi)  ;
robs =    r  ;
labobs = lab + log(Lss) ;
rrate = r - pinfc(+1); 
spreadobsi = exp(rbi(+1))*exp(pinfc(+1)) - exp(robs);
spreadobsc = exp(rbc(+1))*exp(pinfc(+1)) - exp(robs);
//dequity = exp(ga + alfac/(1-alfai)*gv) * ( ( (qcss*scss/lrcss)/(qcss*scss/lrcss+qiss*siss/lriss) ) *(nc-nc(-1)) 
//+ (1- ( (qcss*scss/lrcss)/(qcss*scss/lrcss+qiss*siss/lriss) ) )*(ni-ni(-1)) + z + alfac/(1-alfai)*v );
dtfp = (1-invshare)*(z + (z_l-z_l(-1))) + invshare* (v + (v_l-v_l(-1)))+ (tfp - tfp(-1))  + me_tfpc;

// spread in the data = average spread in the model
// using qx*sx /nx = lrxdata ==> nx = qx*sx *lrxdata^-1
//spread = ( rbc(+1)+log(rbcss)-robs )* (( exp(qc+log(qcss))*exp(sc+log(scss)) )/exp(nc+log(qcss*scss*lrcdata^-1)) ) 
//        + ( rbi(+1)+log(rbiss)-robs )* (( exp(qi+log(qiss))*exp(si+log(siss)) )/exp(ni+log(qiss*siss*lridata^-1)) ) + esp1;
//dcredit = ( qc  - qc(-1) + sc - sc(-1) + z + alfac/(1-alfai)*v + cga + alfac/(1-alfai)*cgv )* ( ( exp(qc(-1)+log(qcss))*exp(sc(-1)+log(scss)) )/ (  exp(qc(-1)+log(qcss))*exp(sc(-1)+log(scss)) + exp(qi(-1)+log(qiss))*exp(si(-1)+log(siss)) ) )
//         + ( qi  - qi(-1) + si - si(-1) + z + alfac/(1-alfai)*v + cga + alfac/(1-alfai)*cgv )* ( ( exp(qi(-1)+log(qiss))*exp(si(-1)+log(siss)) )/ (  exp(qc(-1)+log(qcss))*exp(sc(-1)+log(scss)) + exp(qi(-1)+log(qiss))*exp(si(-1)+log(siss)) ) ) + esp2;

end; 

//steady;

shocks;

var ez;
stderr 0.1611;
var ev;
stderr 1.8366;
var eb;
stderr 1.3419;
var eg;
stderr 0.5019;
var eqs;
stderr 0;
var em;
stderr 0.1175;
var epinfc;
stderr 0.5929;
var epinfi;
stderr 0.2150;
var ew;
stderr 0.3618;
var espi;
stderr 0;

var enc;
stderr 0.0;
var eni;
stderr 0.0;
var ezl;
stderr 0.0;
var evl;
stderr 0.0;
var epki;
stderr 0.0;
var epkc;
stderr 0.0;
var epkc1;
stderr 0.0; 
var epkc2;
stderr 0.0; 
var epkc3;
stderr 0.0;
var epkc4; 
stderr 0.0;
var epkc8; 
stderr 0.0;
var epki1; 
stderr 0.0;
var epki2; 
stderr 0.0;
var epki3; 
stderr 0.0;
var epki4;
stderr 0.0;
var epki8;
stderr 0.0;
var ez1;
stderr 0.1;

var ez4;
stderr 0.1095;
var ez8;
stderr 0.1917;
var ev1;
stderr 0.0;
var ev4;
stderr 0.2352;
var ev8;
stderr 0.7066;
var ethetab;
stderr 0.0;
var etfp;
stderr 0.0;
var etfp4;
stderr 0.0/sqrt(2);
var etfp8;
stderr 0.0/sqrt(2);
var ezl4;
stderr 0.0/sqrt(2);
var ezl8;
stderr 0.0/sqrt(2);
var evl4;
stderr 0.0/sqrt(2);
var evl8;
stderr 0.0/sqrt(2);
var emetfpc;
stderr 0.0;
end;



estimated_params;

// PARAM NAME, INITVAL, LB, UB, PRIOR_SHAPE, PRIOR_P1, PRIOR_P2, PRIOR_P3, PRIOR_P4, JSCALE
// PRIOR_SHAPE: BETA_PDF, GAMMA_PDF, NORMAL_PDF, INV_GAMMA_PDF



stderr ez,1.025,0.01,100,INV_GAMMA_PDF,0.5,2.0;
stderr ev,1.9123,0.001,100,INV_GAMMA_PDF,0.5,2.0;
stderr eb,0.9918,0.025,150,INV_GAMMA_PDF,0.1,2.0;
stderr eg,0.9219,0.01,100,INV_GAMMA_PDF,0.5,2.0;
//stderr eqs,4.27,0.01,120,INV_GAMMA_PDF,0.5,2.0;
stderr em,0.1374,0.01,80,INV_GAMMA_PDF,0.1,2.0;
//stderr espi,0.20,0.01,80,INV_GAMMA_PDF,0.1,2.0;
stderr epinfc,0.046,0.01,80,INV_GAMMA_PDF,0.1,2.0;
stderr epinfi,0.0462,0.01,80,INV_GAMMA_PDF,0.1,2.0;
stderr ew,0.0459,0.01,20,INV_GAMMA_PDF,0.1,2.0;
//stderr enc,0.260,0.01,100,INV_GAMMA_PDF,0.1,2.0;
//stderr eni,0.19,0.01,100,INV_GAMMA_PDF,0.1,2.0;
stderr epkc,2.4268,0.01,100,INV_GAMMA_PDF,0.1,2.0;
stderr epki,4.7757,0.01,100,INV_GAMMA_PDF,0.1,2.0;
//stderr ethetab,0.001,0.0001,10,INV_GAMMA_PDF,0.001,2.0;
stderr ezl,0.3874,0.01,40,INV_GAMMA_PDF,0.5,2.0;
stderr evl,0.3376,0.01,40,INV_GAMMA_PDF,0.5,2.0;
//stderr esp1, 0.05,INV_GAMMA_PDF,0.05,0.01;
//stderr esp2, 0.05,INV_GAMMA_PDF,0.05,0.01;
//stderr epkc1,0.25,0.01,100,INV_GAMMA_PDF,0.1/sqrt(5),2.0;
//stderr epkc2,0.05,0.01,100,INV_GAMMA_PDF,0.1/sqrt(5),2.0;
//stderr epkc3,0.15,0.01,100,INV_GAMMA_PDF,0.1/sqrt(5),2.0;
//stderr epkc4,0.6,0.01,100,INV_GAMMA_PDF,0.1/sqrt(2),2.0;
//stderr epkc8,0.95,0.01,100,INV_GAMMA_PDF,0.1/sqrt(2),2.0;
//stderr epki1,0.05,0.01,100,INV_GAMMA_PDF,0.1/sqrt(5),2.0;
//stderr epki2,0.05,0.01,100,INV_GAMMA_PDF,0.1/sqrt(5),2.0;
//stderr epki3,0.05,0.01,100,INV_GAMMA_PDF,0.1/sqrt(5),2.0;
//stderr epki4,1.03,0.01,100,INV_GAMMA_PDF,0.1/sqrt(2),2.0;
//stderr epki8,0.85,0.01,100,INV_GAMMA_PDF,0.1/sqrt(2),2.0;


//stderr ez1,0.1215,0.01,100,INV_GAMMA_PDF,0.5/sqrt(1),2.0;
stderr ez4,0.1215,0.01,100,INV_GAMMA_PDF,0.5/sqrt(2),2.0;
stderr ez8,0.2287,0.01,100,INV_GAMMA_PDF,0.5/sqrt(2),2.0;
//stderr ev1,0.05,0.01,100,INV_GAMMA_PDF,0.5/sqrt(1),2.0;
stderr ev4,0.1591,0.01,100,INV_GAMMA_PDF,0.5/sqrt(2),2.0;
stderr ev8,1.6733,0.01,100,INV_GAMMA_PDF,0.5/sqrt(2),2.0;

//stderr etfp4,0.13,0.01,100,INV_GAMMA_PDF,0.5/sqrt(2),2.0;
//stderr etfp8,0.52,0.01,100,INV_GAMMA_PDF,0.5/sqrt(2),2.0;


crhoz,0.2673 ,.01,.9999,BETA_PDF,0.4,0.20;
crhov,.0259 ,.01,.9999,BETA_PDF,0.4,0.20;
crhob,.9607,.01,.9999,BETA_PDF,0.6,0.20;
crhog,.9867,.01,.9999,BETA_PDF,0.6,0.20;
//crhoqs,.97,.01,.9999,BETA_PDF,0.6,0.20;
//crhoms,.38,.0000000001,.9999,BETA_PDF,0.4,0.20;
//crhospi,.95,.0000000001,.9999,BETA_PDF,0.4,0.20;
crhopinfc,.6668,.01,.9999,BETA_PDF,0.6,0.20;
crhopinfi,.6677,.01,.9999,BETA_PDF,0.6,0.20;
crhow,.671,.001,.9999,BETA_PDF,0.6,0.20;
//crhosnc,0.82 ,.01,.9999,BETA_PDF,0.5,0.20;
//crhosni,0.55 ,.01,.9999,BETA_PDF,0.5,0.20;
crhospkc,0.2314 ,.01,.9999,BETA_PDF,0.6,0.20;
crhospki,0.7696 ,.01,.9999,BETA_PDF,0.6,0.20;
crhozl,.9966 ,.01,.9999,BETA_PDF,0.6,0.20;
crhovl,.9751 ,.01,.9999,BETA_PDF,0.6,0.20;
//crhothetab,0.7 ,.01,.9999,BETA_PDF,0.5,0.20;
//crhotfp,0.75 ,.01,.9999,BETA_PDF,0.6,0.20;

//stderr emetfpc,0.1044,0.01,100,GAMMA_PDF,0.1, 0.5; 

//cmapc,.28,0.01,.9999,BETA_PDF,0.5,0.2;
//cmapi,.60,0.01,.9999,BETA_PDF,0.5,0.2;
//cmaw,.95,0.01,.9999,BETA_PDF,0.5,0.2;


csadjcost,3.514,0.2,50,GAMMA_PDF,4,1.0;
chabb,0.7650,0.001,0.99,BETA_PDF,0.5,0.1;
//cprobw,0.74,0.01,0.99,BETA_PDF,0.66,0.1;
csigl,0.3441,0.025,10,GAMMA_PDF,2,0.75;
//cprobpc,0.82,0.01,0.99,BETA_PDF,0.66,0.10;
//cprobpi,0.726,0.01,0.99,BETA_PDF,0.66,0.10;
//cindw,0.14,0.000000001,0.99,BETA_PDF,0.5,0.15;
//cindpc,0.13,0.000000001,0.99,BETA_PDF,0.5,0.15;
//cindpi,0.34,0.000000001,0.99,BETA_PDF,0.5,0.15;
czcapi,6.4769,0.000001,50,GAMMA_PDF,5.0,1.0;
czcapc,6.54,0.000001,50,GAMMA_PDF,5.0,1.0;
//czcap,12.2,0.000001,50,GAMMA_PDF,5.0,1.0;
//czcap,0.8,0.000001,1,BETA_PDF,0.5,0.1;

crpi,1.1759,1.0,5,NORMAL_PDF,1.7,0.30;
crr,0.817,0.1,0.99,BETA_PDF,0.60,0.20;
//cry,0.09,0.0000000001,0.8,NORMAL_PDF,0.125,0.05;
crdy,0.1222,0.000001,0.9,NORMAL_PDF,0.125,0.05;

//cminrho,0.17, 0.00000000001,0.999999,BETA_PDF,0.5,0.20;



end;






//varobs dc dinve dy labobs dw  pinfobsc pinfobsi robs ; //spread dcredit;
varobs dc dinve dy labobs dw  pinfobsc pinfobsi robs dtfp; //dequity spread dcredit;
// For Estimation
estimation(datafile= USdata_10_rev, mode_compute=0, mode_file =  twosb_v10b_mobileKNoFItfpNoMe_mode,  prior_trunc=1e-50, first_obs=1, presample=1, lik_init=1, prefilter=0, mh_replic=0,   mh_nblocks=2,mh_jscale=0.25,mh_drop=0.2, nodiagnostic) dc;
//estimation(datafile= USdata_6, mode_compute=4, mode_file = twosb2v6_mode, mode_check, prior_trunc=1e-50, first_obs=1, presample=1, lik_init=1, prefilter=0, mh_replic=0, mh_nblocks=2,mh_jscale=0.20,mh_drop=0.2, nodiagnostic) dc;
//estimation(datafile= USdata_6, mode_compute=0, mode_file = twosb2v6_mode, prior_trunc=1e-50, first_obs=1, presample=1, lik_init=1, prefilter=0, mh_replic=50000, mh_nblocks=2,mh_jscale=0.26,mh_drop=0.2) dc;
//load_mh_file,

// MDD after 50000:
// MDD 



// For IRFs
//estimation(datafile= USdata_6, mode_compute=0, mode_file = twosb2v6_mode, prior_trunc=1e-50, first_obs=1, presample=1, lik_init=1, prefilter=0, mh_replic=0, mh_nblocks=2,mh_jscale=0.25,mh_drop=0.2) dc;
//stoch_simul(irf = 40) dy dc dinve dw pinfobsc pinfobsi robs labobs spreadobsi spreadobsc;

//mode_file =  twosb_v10b_mobileKNoFItfpNoMe_mode, 

//cprobw = 0.001;
//cprobpc= 0.74;
//cprobpi= 0.74;
//cindw=0.5;
//cindpc=0.5;
//cindpi=0.5;
//csigl = 1.63;
//crhoz = 0.9;
//chabb =0.76;

//crpi =2.9;

//czcapc =0.1;
//czcapi =0.1;
//csadjcost = 0.0000000001;

stoch_simul( irf=20) dc dinve invei invec labobs labi labc kpi kpc ;//;dy dinve dc labobs labi labc invei invec spreadobsc spreadobsi qc qi dpi dequity kpc kpi;